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The process of liquid silicon infiltration is investigated for channels with radii from 0.25 to 0.75 [mm] drilled 
in compact carbon preforms. The advantage of this setup is that the study of the phenomenon results to be 
simplified. For comparison purposes, attempts are made in order to work out a framework for evaluating the 
accuracy of simulations. The approach relies on dimensionless numbers involving the properties of the surface 
reaction. It turns out that complex hydrodynamic behavior derived from second Newton law can be made 
consistent with Lattice-Boltzmann simulations. The experiments give clear evidence that the growth of silicon 
carbide proceeds in two different stages and basic mechanisms are highlighted. Lattice-Boltzmann simulations 
prove to be an effective tool for the description of the growing phase. Namely, essential experimental constraints 
can be implemented. As a result, the existing models are useful to gain more insight on the process of reactive 
infiltration into porous media in the first stage of penetration, i.e. up to pore closure because of surface growth. 

A way allowing to implement the resistance from chemical reaction in Darcy law is also proposed. 

Keywords: liquid silicon infiltration, millimeter-sized channels, reaction-formed SiC moiphology, Lattice-Boltzmann simula¬ 
tions for surface growth 


I. INTRODUCTION 

Reactive molten silicon (Si) infiltration into carbon (C) pre¬ 
forms is a process widely employed for the manufacturing of 
ceramic components designed for advanced applications |[lj- 
Ib). The liquid is introduced into the porous structure by capil¬ 
lary forces. This is possible because Si reacts with C to form 
silicon carbide (SiC) GHm. Indeed, Si wets SiC but it does 
not wet C. Another important property of this reaction is that 
the SiC formation leads to surface growth. The thickening of 
the surface behind the invading front results in a retardation 
of the infiltration process and it can also stop the fluid mo¬ 
tion when pore closure occurs. On the other hand, the SiC 
phase is a desired property since providing hardness to the 
material. Furthermore, the surface reaction influences the hy¬ 
drodynamic behavior of the melt significantly, since the infil¬ 
tration dynamics does not follow the usual Washburn behavior 
but it displays a linear dependence on time HEHSl- 

In this study, the experiments of infiltration are performed 
with channels of millimeter width made in dense carbon pre¬ 
forms. This is done with the purpose to gain insight into the 
process of surface growth inside porous media at the micron 
scale. The advantage of the ideal system of a single channel is 
to retain the basic mechanisms underlying the surface reaction 
without the influence of the porous structure. The experiments 
clearly show that SiC formation is determined by diffusion 
first through the liquid and then through the solid SiC layer, 
much slower. Of course, the complexity of the phenomenon 
goes beyond the capabilities of the models for the simulations 
used here (see Sec. m- In order to elucidate this aspect, the 
role of the prior C dissolution and pockets formation is dis¬ 
cussed in detail. Nevertheless, it turns out that the models for 
the simulations proposed in this work can treat the process of 
surface growth leading to pore closure. This is important for 
the optimization of Si-infiltrated ceramic composites and their 
manufacturing. 


The process of infiltration is simulated using the Lattice- 
Boltzmann (LB) method. This approach is based on the dis¬ 
cretization of the velocity space, besides the usual discretiza¬ 
tion of the physical space. The hydrodynamic behavior can be 
reproduced in the incompressible limit at small Mach num¬ 
bers. An increasing variety of problems can be worked out 
exhaustively by this method lEHai. To the end of compari¬ 
son with experiments, we propose a framework for assessing 
the accuracy of simulations. The method is based on dimen¬ 
sionless numbers taking into account the effects of the surface 
reaction. Interestingly, some properties of experimental sys¬ 
tems can be addressed after a suitable calibration. Notably, 
real data for surface growth can be imposed. The outcome 
of simulations proves to be in line with experimental observa¬ 
tions for SiC formation. In general, this work allows to make 
more precise the predictive power of LB simulations for re¬ 
active infiltration thanks to ad hoc experiments. These results 
are valuable in order to sharpen the design and interpretation 
of further research aiming at a direct comparison with experi¬ 
mental work. 


II. MATHEMATICAL BACKGROUND 

Figure [T] shows capillary infiltration into a channel. In 3D 
this phenomenon is described by the differential equation ll22l 

Ea 

[pgiL-z)+ piz]z + 7tr^ (pi - pg) (z)^ = 

2;rr7cos0 — 8;r[pg(L —z)+piz]z . (1) 

The terms on the left hand side come from the inertial forces. 
The first term on the right hand side is due to capillary forces 
and the second one accounts for viscous forces, z is the center- 
line position of the invading front, r and L are the radius and 
the length of the channel, p and p indicate the density and 
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Figure 1: Channel infiltrated by a liquid. The liquid is represented 
in red; in blue the vapor phase. The solid phase is shown in dark. 
For LB simulations, the boundaries of the simulation domain are pe¬ 
riodic, reproducing a flat interface characteristic of a large reservoir. 

the dynamic viscosity; the subscripts 1 and g stand for liquid 
and gas. 7 is the surface tension and 9 the contact angle. The 
behavior described by the well-known Washburn equation is 
obtained in the limit Pg/pi —> 0 and Pg/pi —f 0 . In so doing 
Eq.[T]b ecomes 

(z)2 ^ 27COS0 1 8p . 

Z —- 1 -yZ ■ (2,) 

z pr z pr^ 

A solution to this differential equation is given by ll^ 


exhaustively and efficiently ini. More generally, the LB 
method has an advantage for hydrodynamic studies involv¬ 
ing complex boundaries. In this work we use a multiphase 
model in order to describe the process of infiltration ll25H29l . 
With multicomponent models it is possible to reproduce a lin¬ 
ear behavior for the penetration dynamics ll 22 l[^ . but this is 
not necessary since the weight is the dominant force. In any 
case, it should be kept in mind that multiphase models over¬ 
estimate the experimental density ratio between the gas and 
liquid phases. Previous works addressed the case of linear in¬ 
filtration for single channels IIbTI . 

Specifically, for a d2q9 lattice the dynamics is obtained by 
iterating the LB equation 

fi{r + e,Af,f + Af)=/,■ (r,f) - ^ (r,f(p, m)] 

/ = 0,...,8 . 


z = 


/ yrcos 9 
111 


f+z, 


0 ’ 


(3) 


where zo is the initial position. In the following we will con¬ 
centrate ourselves on this regime and so we will drop the lower 
scripts 1 for liquid. It is worth noticing that this solution suf¬ 
fers from two approximations: first, the contact angle is not 
constant; second, inertia effects are not accounted for. 

Lor channels with radii of millimeter size the weight of the 
liquid contributes to drive the infiltration. This effect is taken 
into account by adding to the right hand side of Lq. |^the ac¬ 
celeration due to gravity g = 9.8 [m/s^]. In our description we 
introduce also the force coming from the column of liquid in 
the reservoir above the channel. Let us assume that the reser¬ 
voir has sides of length Wo, Wi and W 2 with r <C Wb,Wi,W 2 ; 
Wo measures the length along the direction of infiltration z. 
The weight forces deriving from the column of liquid above 
the opening of the channel are described by 


WoWiW2-zr^7: 

W 1 W 2 


r'^npg . 


The fraction gives the height of the column of liquid above 
the capillary. In the course of infiltration, the numerator deter¬ 
mines the new volume of liquid in the reservoir. By dividing it 
by the base area of the reservoir, the new height of the column 
of liquid is obtained. Multiplication by the other terms gives 
the weight of the column of liquid. The above term is added 
to the right hand side of Lq. The differential equation that 
we want to solve is then given by 


(z)^ / 27 COS 0 \ 1 8p . / 

z —-h (-h Wog 1- 2 ^ 3“ ( 1 ~ 

z \ pr J z pr-^ V 

Lor the setup of Lig.[2the last term vanishes. 




W1W2 



This equation is based on the BGK approximation lf3^ . The 
lattice spacing and the time step Af are both assumed to be 1, 
without loosing generality. The vector r = {x,y) points to a 
lattice site. The allowed velocities are indicated by e, and they 
are defined as follows: 


e,- = 


( 0 , 0 ) / = 0 

(cos[(/— l);r/2],sin[(i— l);r/2]) / = 1 —4 

•\/2(cos[(2/ —9);r/4],sin[(2/ —9);r/4]) / = 5 —8 . 


fi is the distribution function for the particles moving with 
velocity e,. In the collision term, the relaxation time is set 
to T = 1 [ts]. The relaxation time determines the kinematic 
viscosity v = (t — l/2)/3. The dynamic viscosity is given by 
p — pv, where p is the density. The equilibrium distribution 
functions read 


fPiP^n) 


Wip{r,t) 


1 + 3ei-u + - 



i = (5) 


The weighting factors are wq =4/9, w,- = 1/9 for i = 1,2,3,4 
and Wi = 1/36 for i = 5, 6 ,7 ,8 ll20l . The local density and 
velocity are computed as follows: 

P{r,t) = Y^fi{r,t) 

1=0 

1 ^ 

Pv P) i=0 


The no-slip boundary condition is implemented by using the 
well-known bounce-back rule. Lluid-fluid interactions give 
rise to interfaces and determine the surface tension 7 . The 
cohesive forces are introduced by means of the function ll25l 

|26l 


III. LB MODELS 



LB simulations are particularly suited for our problem uni- 
En. The reason is that interface phenomena can be treated 


The parameter Gc controls the interaction strength. The func¬ 
tion y/ is defined as (/(r) = ij/o and po are free 
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parameters. The pressure is given by na 

The adhesive forces between the fluid and solid phases deter¬ 
mine the contact angle. They are computed with 


8 

-F’ads(T’) = -Gi^dsWir) Y, + e;Af)e,' . 

r=l 

The parameter Gads sets the interaction strength. The function 
\j/ is the same as for cohesive forces. The function s{x,y) takes 
on the value 1 if the site is occupied by the solid phase, oth¬ 
erwise it is equal to 0. The gravitational force is introduced 
by imposing the constant acceleration g. Finally, the contri¬ 
bution of all the forces is taken into account by using in Eq. 
the velocity u -f t(Fc -f Fads + pg)/P- 

Solute transport takes place on a d2q4 lattice ll34l . with¬ 
out diagonal velocities. As for fluid flow, the evolution of the 
distribution functions is given by Il34l - l38l 

gi{r + eiAt,t+At)=gi{r,t)-^[gi{r,t)-g^‘'{C,u)] 

Ts = 1 [ts] is the relaxation time for solute transport. The dif¬ 
fusion coefficient is given by D = (2 Ts — 1)/4 ll34l . The equi¬ 
librium distribution functions become ll34l 

g^'^{C,u) = ^C[l+2ei-u\ /=l,---,4. (6) 

u is the local solvent velocity; the solute concentration is de¬ 
fined as C = Y.iSi- In order to limit diffusion at the meniscus 
of the invading front we create an interface by using the func¬ 
tion 
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-F’s(t’) =-Gs^(r’)^(p(r + e,Af)e; . (7) 

r=l 

The parameter Gg is used in order to adjust the intensity. The 
function (p is defined as (p{r) = (poe^G/C{r).^ ^ 
arbitrary constants. In analogy with fluid flow, we add to the 
velocity in Eq.|^the term rFg/C. The surface reaction at the 
boundaries is described by the macroscopic equation Il34l438ll 

kr is the reaction-rate constant and Cg is the saturated concen¬ 
tration. In the LB formalism, this condition takes the form 



Figure 2: Sessile droplet for the determination of the contact angle 
6 using LB simulations. The white straight line is the tangent at 
the contact line. The green curve represents the circle fitting to the 
interface. 


H=50 [lu] 



Figure 3: LB simulations for capillary infiltration driven only by ad¬ 
hesive forces. Position of the invading front as time goes on. Points 
for simulation results; solid lines for the theoretical predictions of 
Eq .|^in 2D. L is the length of the capillary and H its height. 


on the solid nodes belonging to the solid-liquid interface is 
obtained by iterating the equation ll37l[3^ 

b{rp+At) =b{r,t)+k^[C — Cf) . 

The initial mass is bo. When the mass b exceeds the thresh¬ 
old value ^maxj the surface grows. This means that one of 
the direct neighbors in the fluid phase becomes a solid node. 
This node is chosen with equal probability among the avail¬ 
able possibilities. The mass of the solid node that triggered 
surface growth is set to zero. 

Hereafter, the output of the simulations will be expressed in 
model units. The symbols for the basic quantities of length, 
mass and time are lu, mu and ts, respectively. In this notation, 
ts stands for timesteps. In the sequel, we will try to com¬ 
pare simulations with experiments by means of dimensionless 
numbers like the capillary number, the Damkohler number 
IIJTII and other numbers involving the surface reaction. This 
procedure guarantees the equivalence of the conditions deter¬ 
mining the hydrodynamic behavior llWll . 


IV. ANALYSIS 


_ +^rCs 

kr+l/2 


and 


Si - -Sx{i) + 2 ^ ■ 


i is the index of the distribution functions that remain unde¬ 
termined after the streaming step; x{i) index for the 

opposite velocity, leaving the fluid phase. Mass deposits on 
the solid surface. Up to a constant, the mass that accumulates 


The Bond number compares gravity and surface tension: 
Bo = i^pgjj. The capillary number compares the effects of 
viscosity and surface tension: Ca = v/f/y, where v is an esti¬ 
mate of the characteristic velocity. In order to compare gravity 
and viscosity, one has to consider the number 


Qx 


Bo 


pp-g 


( 8 ) 


Ca vp 
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Figure 4: LB simulations for infiltration driven only by capillary 
forces. Velocity of the invading front as a function of time. The 
points represent the results from simulations. The solid line is ob¬ 
tained by solving Eq.|^in the 2D case. The dashed lines correspond 
to Eq.[^in 2D. L denotes the length of the channel and H the height. 


L=2500 [lu] 



effects between the reactivity and gravitational forces are thus 
determined by the number 


Higher values of Q 2 mean that reactivity becomes faster 
than mechanics, i.e. the hydrodynamics mainly determined by 
gravitational forces. For completeness, Q^, = prk/jj. compares 
the reaction with the viscosity and Q 4 — prvk/j compares 
the reaction and surface tension. The dimensionless numbers 
are compared in Tab. |^for a typical capillary embedded in a 
porous medium and a millimeter-sized capillary. For r = 10 
[/dm] it is evident that capillary forces drive the infiltration. In 
general, it is assumed that high values for the capillary num¬ 
ber satisfy Ca > 0.01. Qi is high but the viscous resistance 
is in response to capillary forces since Ca^ * is much higher. 
The same conclusions can be drawn by analyzing the other 
dimensionless numbers. In passing to r = 0.5 [mm], gravity 
gains three orders of magnitude with respect to capillarity. As 
a test, if we solve Eq.|^ we find that the term of weight forces 
is 15 times larger than that due to surface tension already after 
2 [sec]. Furthermore, the accordance with Eq.[9]is good. An¬ 
other dimensionless number that is worth mentioning is the 
Damkohler number. It is defined as Da — kH/D, where H 
is a characteristic distance of the systems. Eor example, in a 
purely diffusive regime, this number determines the structure 
of the growing surface; small values for compact structures 
and large values for dendritic morphologies E]. 

Eor a shrinking radius the infiltration kinetics becomes 


Figure 5: Front displacement in the course of time for infiltration 
with a constant acceleration and inert boundaries. Points represent 
the results of LB simulations. For the data in green, the condition 
gH^ = 2-10^^ [lu^/ts^] holds, while the data in yellow satisfy gH^ = 
10^^ [lu^/ts^]. The red solid lines are solutions to Eq.j^in 2D. The 
dashed lines correspond to the predictions of Eq.[T^in 2D. 

We assume that this dimensionless number characterizes our 
systems during infiltration since the gravitational forces find 
the resistance of viscous forces. This means that the resis¬ 
tance of viscosity has a stronger effect than the driving force 
due to adhesive forces. The results of direct simulations with 
Eq. [^indicate that the penetration depth with gravity is always 
smaller than in the case of free fall. As a consequence, if grav¬ 
itational forces are dominant, the typical time scale is set by 
tg = vjg. Note that the velocity v takes into account also the 
effects of viscosity. Eurthermore, under these conditions the 
centerline position of the invading front is given by 


where z{t =0) = zo- This equation is obtained by assuming 
that in Eq.|^gravity and viscosity are the dominant forces, as 
for Eq. 1^ Eurthermore, here reactivity is implemented by im¬ 
posing that the initial radius rq shrinks according to the condi¬ 
tion r = ro — kt, where rq is the initial value ED. The maximal 
penetration depth is attained at the time t — ro/k, leading to 

" 24pk ■ 

A simple calculation indicates that without reaction the infil¬ 
tration depth would be 3zmax during the same time lag with 
zo = 0[lu]. 

If we take into account the weight of the liquid in the reser¬ 
voir as in Eq. for the initial condition of LB simulations 
Eq. [^becomes 


z = 


Pr^8 

8p 


t+ZO ■ 


(9) 


z = 




pr^Lg 

4p 


t+z 


2 

0 ’ 


( 10 ) 


Numerical solutions to Eq.j^with gravity for oil 1401 return re¬ 
sults in good agreement with the above equation. Erom Eq. 
we obtain Qi = 8 for any g realizing the regime of interest. 

The time scale associated with the reactivity is given by 
4 = f/k, where k is the reaction-rate constant. Roughly speak¬ 
ing, this is the time when pore closure occurs. The relative 


with Wo — L the length of the channel. The length of the sim¬ 
ulation domain is given by Nx = 2L. Eor a shrinking radius we 
obtain 

^ = ( 11 ) 
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H=25 [lu] b =810 [mu] 
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H=50[lu] b =610-'[mu] 

max 


0 0.2 0.4 0.6 0.8 I 1.2 1.4 1.6 1.8 2 

time [ts] X lo' 

H=50[lu] b =610“^ [mu] 
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□ g=a 
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H=100[lu] b =6-10"^ [mu] 

max 
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time [ts] X 10 = 

H=100[lu] b =6-10'^ [mu] 
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. 
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Figure 6: Results of simulations in the presence of surface reaction. The length of the capillary is L = 2500 [lu]. The reaction-rate constant is 
chosen so that k = {l/2)H /Ni, where M = 200000 [ts] is the number of timesteps. The first value for the acceleration is given by a = 10^^///^, 
as in Fig.|^ Left: Front displacement as time evolves. Points are associated with simulation results. The red solid lines correspond to the 
theoretical predictions of Eq.[^in 2D. The blue lines represent the predictions of Eq.[T^in 2D. Right: Maximal width of the growing surface 
in the course of time. Points are used for results from simulations. 
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Table I: Characteristic numbers for capillaries of different size infiltrated with molten Si (miiD. Qi compares the relative importance of 
gravity and viscosity, Q 2 of reactivity and gravity, of reactivity and viscosity, of reactivity and surface tension; Bo compares gravity and 
viscosity (Bond number) and Ca compares viscosity and surface tension (capillary number). The reaction-rate constant is set to A: = 4 • 10^* 
[m/s] 1411. T he infiltration velocity for the narrow channel is approximated by v = 1 [/xm/s] 11214141 . For r = 0.5 [mm], the velocity is estimated 
using Eq.M The surface tension of Si is y = 0.86 [N/m], the dynamic viscosity /t = 0.94-10^^ [N-s/m^] and the density p = 2.53 • 10^ [kg/m^] 
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The maximal penetration depth is in this case 


^max — 




Yl\xk 


+ z 


During the same time interval, without reaction the infiltra¬ 
tion depth would be v/Szmax by assuming that zo = 0 [lu]. The 
above results can be useful in order to obtain estimates for 
the characterization of experiments. For the sake of complete- 
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Q2=0.02 



Figure 7: Equivalence of systems based on the relative strength be¬ 
tween reaction and gravity, determined by the characteristic number 
Ql = vk/{rg). The results are plotted for the 2D case according to the 
representation suggested by Eq.[T^ In this specific case, the quantity 

I is given by I = gtr^/k. For the acceleration we use a = 10^^///^; 
Ni = 200000 [ts] is the total number of timesteps 


ness, if capillary forces remain significant it is obtained 


z = 


['■q - ('■o - ktf] + VI - ('■o - ktf] +zl . 


Uiik^ ^ J Afik 

The maximal penetration depth becomes 


( 12 ) 


Zmax — 


I pLgrl yrl cos 0 


+ Z| 


0 ■ 


I2jik 4pk 

Equation[^can be rewritten in dimensionless form as 


(13) 


where t = kt /kq and z = z/l with I = ^J~pLgr^J{pik) ■ Instead, 
if capillary forces are the dominant ones, Eq.[T^can be recast 
into dimensionless form as 


cos 0 


z -Zo = 


[l-(l-t) 2 ] 


where t =kt/ro and z = z/l with I = tq ^yj(jik). 


B. Contact angle 

The equilibrium contact angle 0 is obtained by using the 
method of the sessile droplet. A droplet of liquid is put in 
contact with the surface. The droplet spreads in order to reach 
the equilibrium state consisting of a spherical cap. The tan¬ 
gent at the contact line, where the three phases meet, defines 
the contact angle (see Eig.[^. Eor adhesive forces i^ads, the po¬ 
tential parameters are given by Gads = —380 [lu/ts^], Xj/Q — A 
and po = 200 [mu/lu^]. The outcome for the equilibrium con¬ 
tact angle is 0 = 32°. The experimental value for the contact 
angle is 30° lfT2l[T3]l. 

C. Washburn infiltration 

In experiments based on pure Si and some alloys with 
porous C, the infiltration process exhibits a linear time de¬ 
pendence for the displacement of the invading front OMl. 
The breakdown of the ordinary Washburn behavior of Eq. |^is 
ascribed to the effects due to the reaction between Si and C to 
form SiC lfT2Ul4l . Since gravity is expected to be the domi¬ 
nant force, in this study we use models reproducing the usual 
parabolic time dependence ll25H29l . In our systems, the liquid 
to vapor density ratio is 10, so quite high. Eirst, we thus eval¬ 
uate the implications of this approximation. Simulations for 
the process of infiltration are performed with systems similar 
to that shown in Eig.[T] Eigurej^ shows the dependence of the 
results on the length of the capillary L for a given height H. It 
turns out that, as the height increases, the departure from the 
predictions of Eq.|^is more pronounced. The dependence on 
the length of the channel is less marked. It is important to re¬ 
call that what characterizes the hydrodynamic behavior is the 
velocity. Eigure|^ compares the velocity to the predictions of 
the models. Also for H = IQO [lu] the velocity is almost equal 
to that obtained with a vanishing density for the gas phase 
relatively soon. This implies that the simulated systems re¬ 
produce correctly the desired Reynolds number and capillary 
number of the analytical and numerical solutions. As a result, 
in the sequel simulations of validation will be carried out with 
capillaries of height at most of 100 [lu]. 


D. Capillary infiltration with gravity 


V. SIMULATION RESULTS 
A. Surface tension 

To start with, we determine the surface tension y using 
Laplace law. This law states that, for a liquid bubble immersed 
in its vapor phase, the pressure drop across the interface is pro¬ 
portional to y. In 2D, Laplace law reads AP ^ Pin-Pout = 7/r. 
Eor cohesive forces Fc, the potential parameters are set to 
Gc = —140 [mu-lu/ts^], 1 / 0=4 and po = 200 [mu/lu^]. It 
is found that the surface tension is 7 = 29.91 [lu mu/ts^]. The 
densities for the liquid and the gas are pi = 700 [mu/lu^] and 
Pg = 70 [mu/lu^]. The density ratio is thus Pi/Pg = 10. 


We consider capillary infiltration in the presence of a con¬ 
stant acceleration mimicking gravity. Eigure shows the re¬ 
sults for LB simulations. The two groups of data are charac¬ 
terized by different values for the product gH^. As predicted 
by Eq. we see that within a group of results the behavior 
of the systems is comparable (cf. Eq. 10 1 , with the exception 
for the initial stages of the dynamics. The deviations are more 
marked for increasing width of the channels and for increas¬ 
ing acceleration. Under these conditions the viscous and in¬ 
ertia resistance of the gas phase results to be stronger. The 
main effect is that at the beginning the process is particularly 
slow. As a result, we consider the present systems reliable for 
the investigation of the properties of surface growth. Table [n| 
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//[lu] 

Qi 

Bo 

Ca 

25 

1.9 

5.7-10-2 

2.9- 10-2 

50 

1.8 

5.8-10-2 

3.2- 10-2 

100 

1.6 

5.8-10-2 

0 

1 

to 


Table II: Characteristic numbers for infiltration driven by gravity 
without surface reaction (see Fig. [^. The constant acceleration is 
g = 10^^/H^. H is the width of the channel. Qi compares gravita¬ 
tional and viscous forces. Bo compares gravity and surface tension 
(Bond number), while Ca compares viscosity and surface tension 
(capillary number). The quantities that can vary in the course of time 
are averaged over a time interval centered around Nt/2 = 100000 
[ts]. 


H [lu] Qi Q 2 _ Q 3 _ 04 _Bo_ Ca 

25 0.87 2.9-10^^ 2.5-10^^ 5.4-10^^ 1.9-10^2 2.2-10^2 

50 0.96 1 . 1 - 10-2 1 , 0 . 10-2 2 . 2 - 10-4 2.1-10-2 2.2-10-2 

100 0.93 4.6-10-2 4.3-10-2 1.1-10-2 2.5-10-2 2.6-10-2 


Table III: Characteristic numbers for simulations with constant ac¬ 
celeration and surface reaction (see Fig. |^. H is the initial width 
of the channel. Qi compares gravitational and viscous forces. Q 2 , 
Qt, and Q4 compare the effects of reaction to gravitational, viscous 
and capillary forces. Bo is the Bond number comparing gravity and 
surface tension; Ca is the capillary number comparing viscosity and 
surface tension. The results for a given initial height are averaged 
over the different applied accelerations (see Fig. |^. The quantities 
that can vary as time evolves are averaged inside a window centered 
around the time Ni/2 = 100000 [ts]. 


reports the characteristic numbers of relevance for these set¬ 
tings. As indicated by the simulation results, these values are 
associated with a dominant force given by gravity. 


//[lu] 

Qi 

02 

23 

24 

Bo 

Ca 

25 

1.3 

1.3-10-2 

1.7-10-2 

4.6-10-2 

3.6-10-2 

2.7- 10-2 

50 

1.4 

4.9-10-2 

7.0-10-2 

2.0-10-4 

4.2-10-2 

2.9- 10-2 

100 

1.3 

2.2-10-2 

2.9-10-2 

9.6-10-4 

4.5-10-2 

3.3-10-2 


Table IV: Characteristic numbers for simulations with constant ac¬ 
celeration and reactivity enabled. The reaction-rate constant is fc = 
(1/4)///M. The results thus fulfill the condition r[t = 200000) s:; 
///4. H is the initial width of the channel. Q\ determines the relative 
strength between gravitational and viscous forces. Q2, 23 and Q4 
compare the reaction with gravitational, viscous and capillary forces. 
Bo is the Bond number comparing gravity and surface tension; Ca is 
the capillary number comparing surface tension and viscosity. For 
a given initial height an average is taken over the different applied 
accelerations, chosen as in Fig. The time-varying quantities are 
averaged inside a window centered around the time Nt/2 = 100000 
[ts]. 


Boat shape crucible Sample Silicon Silicon carbide cloth 



Figure 8: Top: Designed setup for experiments. The description 
and function of the indicated parts can be found in main text (see 
Sec. I VI a]). Bottom: Schematic for the sample with drilled channels. 



Figure 9: Micrograph near the inlet for the channel of diameter 1 
[mm] infiltrated with a dwell time of 5 [min]. 


E. Infiltration driven by weight with surface reaction 

We now enable the reaction leading to surface growth. 
The initial conditions for solute transport are Ci = 1.2 -10^^ 
[mu/lu2] in the liquid phase and C 2 = 7.5 -10^2 [mu/lu2] in the 
vapor phase. The saturated concentration is set to Cs = 5 -10^2 
[mu/lu2]. The parameters for the function of Eq.|^ introduc¬ 
ing an interface are given by Gs = —4.875 -10^2 [mu/lu/ts2], 
(po = 1 and Co = 4.9 - 10^2 [mu/lu2] 1^1 [H. The initial 
mass on the solid boundaries is ho = 2 - 10^2 [mu]. To start 
with, the reaction-rate constant is set tok= {l/2)H/Nt, where 
At = 200000 is the total number of iterations. This condition 
alone is not sufficient in order to ensure that pore clogging oc¬ 
curs at time At, since there remains to be fixed the threshold 
value hmax for the cumulated mass. Figurej^shows the results 
for the suitable values of hmax- For two systems the surface 
grows more slowly than expected. This turns out to be related 
to the formation of the interface for solute near the opening 
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Figure 10: Cross-section of infiltrated channel of diameter 0.5 [mm] 
with dwell time of 24 hours. The circle represented with a solid line 
indicates the initial diameter. The process of C dissolution leads to 
an increase of the diameter to 9.2 ± 4 [/xm] on average. 



Figure 11: Formation of pockets in the growth of SiC. The initial 
diameter of the channel is 0.5 [mm] and the experiment proceeded 
with a dwell time of 24 hours. 


of the channel preventing solute from reaching properly the 
solid surface. This shortcoming can be reduced by consid¬ 
ering higher initial values for the solute concentration in the 
liquid phase. For comparative purposes we keep the present 
settings. In Tab. [HI] are listed the results for the characteristic 
numbers. The values for 03 ^nd Qa, are representative for 
systems dominated by reactivity, as it is evident from Fig. [^ 
The reactivity develops a stronger influence as the width of the 
channel doubles with an increase of Qi, Qi and Qa, approxi¬ 
mately by a factor of 4. Gravity remains the dominant force: 
an estimate of the terms due to gravitational and capillary 
forces in the differential equation reveals that the latter is 9.6 
times smaller (see Eq.[^. Comparison with the results without 
reaction. Tab. [^ indicates that the viscosity becomes stronger 
with respect to gravity. Furthermore, surface tension becomes 



SiC grains boundaries 


Figure 12: Schematic representation of the two separate stages at 
the origin of SiC growth with the formation of pockets 1431 . 



2 3 4 

channel length [mm] 


Figure 13: Profile of the SiC layer as resulting from simulations 
after 12.5 [min]. As the height of the channels triples, the width of 
the reaction-formed phase can be assumed to remain comparable, as 
in experiments. Top: upper wall; Bottom: lower wall. H is the initial 
height of the channels, i.e. the initial diameter. The initial length of 
the channels is L = 7 [mm]. 


more relevant with respect to gravity. We thus conclude that 
the theoretical predictions shown in Fig.|^ including also cap¬ 
illary effects, are closer to the simulation results for this rea¬ 
son. Actually, we could estimate that for H = 50 [lu], in the 
absence of reactivity, the term due to surface tension is 15 
times smaller than that introducing a constant acceleration. It 
is also found that the spatial variations of the density for the 
simulations of Fig.j^are comparable to those of Fig. [^ where 
the accordance with the theory is reasonable. This means that 
the treatment of the forcing term does not represent a signif¬ 
icant source of inaccuracy ll42l . In Fig. [^ the infiltration re¬ 
sulting from simulations is still faster presumably because it 
is assumed that the channels shrink uniformly. The fact that at 
the contact line the width of the channel is the initial one is in 
particular relevant in 3D (see Eq. [T]). The results for H — 50 
[lu] clearly show that surface growth can not be disregarded 
when the reduction of the radius is larger than 1 /2. Thus, we 
consider surface reaction as the dominant phenomenon. The 
best precision is obtained for H = 50 [lu], when inertia forces 
are still weak at the beginning. In Tab. |In|we report the char¬ 
acteristic numbers for the condition k = {l /4)H/Nt. In words, 
this means that at time At, the end of the simulation, the initial 
radius is divided by 2. It can be seen that gravity becomes 
more relevant than in the previous case as indicated by 2j, (22 
and Bo. Also viscosity becomes more significant with respect 
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to the reaction. The relative strength among the other forces 
remains almost unchanged. Again, Q 2 , ^3 and 04 increase 
by a factor of 4 as the initial height doubles. In Fig. we 
verify the equivalence of systems based on the characteristic 
number Q^. The agreement between the curves is good. This 
representation confirms that the reactivity and gravity account 
for the dominant effects. As discussed before, the discrepancy 
with the theoretical results can be ascribed to the fact that cap¬ 
illary forces remains significant. It should also be noted that in 
the representation of the results of Fig. the difference turns 
out to be of 3%. 


VI. EXPERIMENTS 

A. Material, setup and method 

Figure shows the setup egineered for the experiments 
of reactive infiltration into millimeter-sized channels. The 
graphite used in order to prepare the preforms was provided 
by Schunk Graphite Technology. The porosity of the sam¬ 
ples is constrained below 6 %, the density is 1.94 [g/cm^] with 
very low ash content. As-received Si is a powder consist¬ 
ing of flakes at most of 500 [fim] with a purity of 99.9%. 
Fine control of the infiltration process is achieved by making 
use of alumina boat-shaped crucibles for the assembly of the 
components (see Fig. [^. High-density alumina is a common 
choice for high-temperature technologies due to its thermal, 
mechanical and chemical resistance. A SiC cloth is prepared 
as support for Si flakes and in order to filtrate Si oxides at 
the time of infiltration. The auxiliary elements are in graphite 
coated with boron nitride so as to prevent molten Si from re¬ 
acting. Undesired effects of residual oxygen are avoided by 
surrounding the sample with activated carbon particulate with 
size less than 1 [mm]. The experiments are performed in a 
horizontal tube furnace at a temperature of 1450° [C] (heat¬ 
ing rate of 3° [C/min]) under argon atmosphere (flow rate of 
100 [cm^/min]). The infiltrated samples are sectioned and pol¬ 
ished for surface analysis with optical microscopy. The chan¬ 
nels are realized in the compact graphite with a drill press or 
with a milling machine for the smallest diameter of 0.5 [mm]. 
Besides compressed air, the cleaning treatment of the chan¬ 
nels with smallest diameter consists also of ethanol ultrasonic 
bath. All the channels have a length of 7 [mm]. The diameter 
of the channels and the dwell time vary in the experiments. 


B. Experimental results 

The first set of experiments include channels of diameters 
of 0.7, 1 and 1.5 [mm] for hold times of 1 hour, 30 [min] 
and 5 [min]. It turns out that the microstructure of all samples 
presents the same morphology and the same rate of growth for 
SiC. It is observed that the width of the reaction-formed SiC 
phase varies between 6 and 25 [/im] (cf. Fig. |^. Inspection 
of the results indicates that the thickness of the SiC layer is 
independent of the capillary diameter. 


In the second series of experiments the diameter of the 
channels is fixed to 0.5 [mm] and the dwell time is increased 
to 24 hours. The motivation for such a dwell time is to verify 
the slow process of SiC formation following an initial growth. 
The obtained channels have a good circularity with the av¬ 
erage diameter given by 516 ± 3 [/im]. After infiltration the 
diameter increases to an average value of 525.5 ±1.5 [/tm], as 
it can be seen in Fig. This change in diameter of the chan¬ 
nels can be associated with the C dissolution process during 
the initial stage ca. Also in this case, the crystals of SiC at¬ 
tain a thickness between 6 and 21 [fim]. The only difference 
that can be identified is that the formed SiC is more compact 
for longer dwell times. It is to be noted that SiC grows up¬ 
ward leaving behind its iterface regions filled with Si, usually 
referred to as pockets. It turns out that during cooling, SiC 
crystals of submicrometric size precipitate. Owing to C dis¬ 
solution, the thickness of the SiC layer can reach 38 [/tm], as 
shown in Fig. 11 

In another set of experiments the diameter of the channels is 
again 0.5 [mm] and the dwell time is set to 72 hours. We point 
out that the SiC layer is not uniform along the wall of the chan¬ 
nels. It is found that near the inlet the thickness is between 20 
and 30 [fim]. By advancing into the channels, the SiC thick¬ 
ness decreases up to reaching a constant value in the range of 
10—15 [/tm]. As before, the diameter of the capillaries after 
infiltration increases to an average value of 524 ± 1 [/tm]. In 
general, the growth of the SiC layer is similar to that obtained 
with a dwell time of 24 hours. The main significant differ¬ 
ence is that the SiC phase becomes more compact. Another 
differing aspect is that the pockets become lower in number 
and more compact. This observation corroborates the mech¬ 
anism of SiC formation proposed in Ref. Ii43l . Namely, SiC 
growth proceeds in two steps. At the beginning, the growth 
is faster because of the nucleation of SiC crystals close to the 
C surface (distance in the range of 10 — 20 [/ 4 m]). Between 
these nuclei. Si dissolves C which is transported through the 
liquid to the growing crystals. This explains the quite regu¬ 
lar distance between pockets formed at SiC-C interface. This 
process is shown in Fig. [T^ The second phase starts when 
the crystals merge. At this point the process is governed by 
the mechanism of Si and C diffusion through SiC, slower than 
in the case of transport through the liquid after C dissolution 
by Si. The result is that the overall reaction rate slows down 
drastically. This is why the difference in thickness of SiC is 
slight for experiments performed during 24 and 72 hours. 


C. Comparison with LB simulations 

Equation points out that in experiments the infiltration 
driven by gravity is a transient process. Also in the exper¬ 
imental work, the column of fluid above the opening of the 
channels can be assumed to be comparable to the length of 
the channels. As a consequence, the surface grows mainly un¬ 
der a regime of diffusion for mass transport. This means that 
an equivalence with the simulations can be established via the 
Damkohler number Da ifJTll . In order to estimate the experi¬ 
mental Da numbers we consider for the characteristic length 
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of the systems the side of an equivalent square. Given the dif¬ 
fusion coefficient D = 5 • 10^® [m^/s] IH] and the reaction- 
rate constant k = A- 10^^ [m/s] of Si, it is found that the Da 
number varies in the range 0.015 — 0.026. It follows that the 
process of Si diffusion is much faster than the reactivity. In the 
LB systems, the length of the channels is set ioL— 1400 [lu], 
corresponding to 7 [mm] (see Sec. VI Aj i. Instead, the width 
of the simulation domain is determined by the aspect ratio of 
the experimental channels, leading to Ny = 100,140,200,300 
[lu]. By reproducing the experimental Da numbers, it is ob¬ 
tained that the reaction-rate constant for the simulations is the 
same for every system and equal to k — 10^^ [lu/ts]. In the 
LB models for the reaction, the surface grows following a lin¬ 
ear relation with time. For this reason we can only simulate 
the process up to the formation of a continuous SiC layer (see 
Fig. 12 1 . The number of timesteps is derived by imposing that 
the surface grows for a distance corresponding to 30 [/rm] dur¬ 
ing this time. It is found that the iterations amount to 6 • 10^ 
[ts]. The simulation domain is 3000 [lu] long. The channels 
start at X = 1000 [lu] {L — 1400 [lu]). The system is filled 
with liquid and solute up to x = 2500 [lu]. The channel is thus 
immersed in the liquid and the system remains at rest. This 
state could have been reached after a previous phase of infil¬ 
tration driven by an external acceleration. In experiments this 
regime is very short (see Eqs.[9|and[T0|i. The motivation of 
our choice is to maintain the same conditions for mass depo¬ 
sition. Indeed, tests show that a contact angle of 90° results 
in the possible appearance of the interface for solute near the 
walls. 


Figure [T3] shows the final profile of all channels. It is inter¬ 
esting to note that, as the height triples, the results remain 
comparable. Furthermore, the target value of 30 [/im] for 
the maximal width of the SiC layer is always neared. These 
properties are in agreement with experiments. For the larger 
channels, the process of surface growth is slower. A possible 
reason is that the walls are more distant from the centerline, 
where the concentration of solute is higher. On the other hand, 
in this case, the SiC layer tends to be uniform because there is 
more solute available in the systems. The simulations fail to 
reproduce a profile decreasing inside the channels. It is found 
that the surface grows from the extremities to the center. As a 
consequence, a factor responsible for this shortcoming is that, 
in simulations, there is an additional source of solute close to 
the right opening of the channels. It turns out that this phe¬ 
nomenon is less marked with the smaller channel. In systems 
for porous structures the channels are quite narrower and not 
blind. These undesired effects should be weaker also by the 
presence of advection induced by the advancing front. 

Additional experiments are carried out in order to single 
out the causes leading to a discrepancy for the SiC profile. 
The channels have a diameter of 1 [mm] with openings of in¬ 
let and outlet, as in simulations. Furthermore, the preforms 
are in glassy carbon since this carbon is characterized by a 
reduced roughness, porosity and reactivity (Alfa Aesar sup¬ 
plier). The experimental conditions are the same as before, 
except for the temperature set to 1475° [C] and the size of the 
Si flakes (below 5 [mm]). In the first experiment, the length 
of the channel is 7.85 [mm] and the holding time is 2 hours. 


At the outlet, the flow is stopped by the formation of SiC crys¬ 
tals. It is found that the SiC growth around the capillary wall 
is uniform, forming a layer 25 — 35 [fim] thick. It turns out 
that the microstructure of the reaction-formed phase exhibits a 
higher roughness. In the second experiment, the length of the 
channel is 5.71 [mm] and the holding time is 4 hours. In this 
case, at the outlet the infiltration stops with the formation of 
a pendant droplet. It is observed that the SiC layer formed of 
is more compact but there is no significant difference in thick¬ 
ness. Again, both phenomena can be ascribed to the slow pro¬ 
cess of C diffusion through the SiC layer, affecting only the 
density. These results for channels with outlet compare better 
with the simulations of Fig. 13 They also prove that the initial 
fast infiltration under the action of gravity has a marginal role 
for surface growth, as assumed in simulations. 


VII. RESISTANCE FROM CHEMICAL REACTION IN 
CAPILLARY INFILTRATION 


Let us indicate by F any force acting on the meniscus, 
i.e. the interface, in Eq. Let us assume that the radius 
shrinks as r(f) = ro — kt. Erom the experiments we also know 
that z = vt ESlIIl. If we replace into Eq. with inertia ne¬ 
glected, it is found that the resistance from the reaction has 
the explicit form 03 


F 


—iTtroycos 9 + 8n7t(v^ + 


kycos9 

4/1 


t . 


(14) 


Here v is the interstitial velocity for a single channel. 

Eluid flow through a porous medium is generally treated 
using the empirical result known as Darcy law. It states that 
the following relation holds 


_ KAP 
p z 


(15) 


K is the permealibity and P the pressure. Here u is the super¬ 
ficial velocity; by dividing it by the porosity, the interstitial 
velocity is obtained (continuity of fluid flow). In the Hagen- 
Poiseuille model, the porous medium is regarded as an assem¬ 
bly of parallel capillary tubes IBl. In this model, for reactive 
flow, the permeability takes the form BTl l44ll 


(PO , <P{t) ,^.2 

K=^Jn-b) = —r{<) 


(po is the initial porosity. The porosity in the course of time is 
given by (p(t) = (po{ro — ktY/r^. In order to integrate Eq. 
we also use the relation 


15 


3kt 

(poro[l-(l-^f)']^’ 


q is the distance allowing to calculate the fluid flow in the 
reservoir. When pore closure occurs z = 3q/<pQ, where <Pq/3 
is the average porosity. 
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If the resistance arising from the chemical reaction of 
Eq. [T^is taken into account, in Eq. [^the pressure becomes 
AE =/TTr^. Darcy law, Eq. 


15 


leads to 


2 _ fty ^ 


i-(i-At)^ 

V ro / 


(16) 


The velocity v is a representative value for the velocity inside 
the capillary tubes. The dependence on the surface tension, 
the contact angle and the viscosity is incorporated into the 
velocity v. The infiltrated distance is determined by the fol¬ 
lowing integral (q derivative with respect to time) 


Z-ZO = 


Jo (p 

(porov‘ 

3k 


-dt 


f 

Jo 


l-(l-^f)3 

^ ro ^ 


-dt 


<po^y 

9k^ 


[i-(i-^0 'P+?^ 


If q'o = 0 we obtain the expected result z = vt and thus Eq. 
is correct. 


VIII. CONCLUSIONS 

In this study, the focus is on millimeter-sized channels in 
order to elucidate relevant aspects of Si infiltration into C pre¬ 
forms. A systematic series of experiments is performed with 
the purpose to detail the process of surface growth. The ex¬ 
perimental work clearly shows that the formation of SiC from 
the reaction between Si and C proceeds through two stages. 
The first one is mainly responsible for the thickening of the 
SiC layer. This stage is characterized by diffusion of reactants 
in the liquid. In the second stage, the diffusion of reactants oc¬ 
curs through the reaction-formed SiC solid phase. This mech¬ 
anism is much slower and has the effect to make the SiC layer 
more uniform and compact. 

The system of a single channel has also the advantage to 
ease the comparison with simulations in order to improve the 


description of the surface reaction MM- Dimensionless 
numbers are introduced with the aim to assess the accuracy 
of simulations involving surface growth. Our analysis points 
out that the characteristic numbers lead to a consistent ap¬ 
proach. This is critical for extensions to porous structures. 
Eurthermore, this work indicates that the models for the re¬ 
action support a proper implementation of the pore-closing 
phenomenon. This means that the important parameter of the 
reaction-rate constant can be regarded as an effective rate of 
growth. Constraints arising from experiments can be imposed 
for the profile of the SiC layer. There appears that in the diffu¬ 
sive regime the structure of the growing phase can be regarded 
as a crude approximation of the microstructure of SiC. It is to 
be noted that infiltration velocities are low for real porous pre¬ 
forms insiiia. 

Einally, this work represents a substantial step toward the 
simulation of reactive infiltration into porous structures. Our 
investigation suggests that a correspondence between infil¬ 
trated distances in simulations and experiments can be estab¬ 
lished. Eurthermore, the comparison can be made more strin¬ 
gent by means of characteristic numbers taking into account 
the surface reaction. The properties of the forming SiC layer 
reproduced by the models are also a positive indication. The 
prospect is to address the involved problems of maximal pen¬ 
etration, SiC formation and residual Si content for the opti¬ 
mization of ceramic materials. 
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